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Abstract 

The computation of the energy spectra of Standard Model particles originated from 
the annihilation/decay of dark matter particles is of primary importance in indirect 
searches of dark matter. We compute how the inclusion of electroweak corrections 
significantly alter such spectra when the mass M of dark matter particles is larger 
than the electroweak scale: soft electroweak gauge bosons are copiously radiated 
opening new channels in the final states which otherwise would be forbidden if 
such corrections are neglected. All stable particles are therefore present in the final 
spectrum, independently of the primary channel of dark matter annihilation/decay. 
Such corrections are model-independent. 
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1 Introduction 

There are overwhelming cosmological and astrophysical evidences that our universe contains a 
sizable amount of Dark Matter (DM), i.e. a component which clusters at small scales. While its 
abundance is known rather well in terms of the critical energy density, ^2dm^ 2 = 0.110±0.005 [1], 
its nature is still a mistery. Various considerations point towards the possibility that DM is 
made of neutral particles. If DM is composed by particles whose mass and interactions are 
dictated by physics in the electroweak energy range, its abundance is likely to be fixed by the 
thermal freeze-out phenomenon within the standard Big-Bang theory. DM particles, if present 
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in thermal abundances in the early universe, annihilate with one another so that a predictable 
number of them remain today. The relic density of these particles comes out to be: 

tt DM h 2 _ 3 x 10- 26 cm 3 /sec 

0.110 (av} ann ' U 

where (cTv) aim is the (thermally-averaged) cross annihilation cross sections. A weak interaction 
strength provides the abundance in the right range. This numerical coincidence represents the 
main reason why it is generically believed that DM is made of weakly-interacting particles with 
a mass in the range (10 2 — 10 4 ) GeV. There are several ways to search for such DM candidates. 
If they are light enough, they might reveal themselves in particle colliders, such as the LHC, 
as missing energy in an event. In that case one knows that the particles live long enough to 
escape the detector, but it will still be unclear whether they are long-lived enough to be the 
DM [2]. Thus complementary experiments are needed. In direct detection experiments, the 
DM particles elastically scatter off of a nucleus in the detector, and a number of experimental 
signatures of the interaction can be detected [3]. In indirect searches DM annihilations or decays 
around the Milky Way can produce Standard Model (SM) particles that decay into e ± ,p, p, 7 
and d , producing an excess in their cosmic ray fluxes. Present observations are approaching 
the sensitivity needed to probe the annihilation cross section suggested by cosmology, eq. (1). 

Furthermore, this topic recently attracted interest because the PAMELA experiment [4] 
observed an unexpected rise with energy of the e + / (e + + e _ ) fraction in cosmic rays, suggesting 
the existence of a new positron component. The sharp rise might suggest that the new com- 
ponent may be visible also in the (e + + e~) spectrum: although the peak hinted by previous 
ATIC data [5] is not confirmed, the FERMI [6] and HESS [7] observations still demonstrate a 
deviation from the naive power-law spectrum, indicating an excess compared to conventional 
background predictions of cosmic ray fluxes at the Earth. While the current excesses might be 
either due to a new astrophysical component, such as a nearby pulsar [8], or to some experi- 
mental problem, it could be produced by DM with a cross section a few orders of magnitude 
larger than in eq. (1), maybe thanks to a Sommerfeld enhancement [9, 10]. 

In any case, it is undeniable that nowadays indirect search of DM is a fundamental topic in 
astroparticle physics, both from the theoretical and experimental point of view. Computing the 
energy spectra of the stable SM particles that are present in cosmic rays and might originate 
from DM annihilation/decay is therefore of primary importance. 

The key point of this paper is to show that electroweak radiative corrections have a siz- 
able impact on the energy spectra of SM particles originated from the annihilation/decay of 
DM particles with mass M somehow larger than the electroweak scale. The reason is in fact 
simple and should be familiar to readers working in collider physics: at energies much higher 
than the weak scale (in our case the mass M of the DM) soft electroweak gauge bosons are 
copiously radiated from highly energetic objects (in our case the initial products of the DM 
annihilation/decay). This emission is enhanced by In M 2 /M^ when collinear divergences are 
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present and In 2 M 2 /M' 2 V when both collinear and infrared divergences are present [11]. These 
logarithmically enhanced terms can be computed in a model-independent way through the well 
known partonic techniques based on the Collinear Approximation (CA). Our work will involve 
generalizing the partonic splitting functions to massive partons, because our 'partons' include 
the W, Z bosons. 

Putting these technical details aside, what is important is that the emission of gauge bosons 
changes significantly some final energy spectra. Indeed, suppose that DM annihilates into a 
pair of leptons. The emitted gauge bosons give hadrons (resulting in a p flux) and mesons 
(giving a significant extra amount of photons via 7r° — > 77). The total energy gets distributed 
among a large number of lower energy particles, thus enhancing the signal in the lower energy 
region (say, (10 — 100) GeV), that is measured by present-day experiments, like PAMELA. 

This paper will be inevitably rather technical and therefore we have decided to defer as 
many as possible technicalities to the various Appendices. To diminish the burden, we present 
qualitative considerations in Section 2 and outline the quantitative computation in Section 3. 
The reader interested in the final results may jump directly to Section 4 where our findings are 
presented. Conclusions are presented in Section 5. In Appendix A we discuss EW evolution 
equations, in Appendix B we derive parton splitting functions for massive partons, and in 
Appendix C we list all splittings among SM particles, including the effect of the top Yukawa 
coupling. 

2 Qualitative discussion 

As mentioned in the Introduction, the presence of DM is probed indirectly by detecting the 
energy spectra of stable particles (p ± ,e ± ,u, r y,d). At a first sight, since electroweak radiative 
corrections are expected to be small — weak interactions are weak, after all — they might 
seem to play no role in the DM indirect searches. At the typical weak scale of (9(100) GeV, 
radiative corrections produce relative effects of (9(0.1)%. For instance, this was the case for 
experiments that took place at the LEP collider. However at energies of the order of the TeV 
scale, like those probed at the LHC, things are different: electroweak radiative corrections can 
reach the (9(30)% level [12] and they grow with energy, eventually calling for a resummation 
of higher order effects [13]. In a nutshell, what happens at energies much higher than the 
weak scale is that soft electroweak gauge bosons are copiously radiated from highly energetic 
objects that undergo a scattering with high invariant mass. This is much the same as photon, 
or gluon, radiation whenever the hard scale is such that the W, Z masses can be safely taken 
to be very small. Important differences with respect to unbroken gauge theories like QCD and 
QCD arise in the case of a (spontaneously) broken theory like the EW sector of the SM. It was 
found [14] that in hard processes with at least two relativistic non abelian charges, effective 
infrared divergencies that are manifest as double log corrections («2 hi 2 M 2 /Myy) appear. They 
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Figure 1: DM annihilation/ decay initially produces a hard positron- electron pair. The spectrum 
of the hard objects is altered by electroweak virtual corrections (green photon line) and real Z 
emission. The Z decays hadronically through a qq pair and produces a great number of much 
softer objects, among which an antiproton and two pions; the latter cascade decay to softer 7s 
and leptons. 

are not present in QED and QCD and this effect has been baptized "Bloch-Nordsieck Theorem 
Violation" [14]. We refer the reader to the relevant literature [14, 15, 16] for details. In the case 
at hand, since the initial DM particles are nonrelativistic, radiation related to the initial legs 
does not produce log-enhanced terms. Therefore, we only need to examine soft EW radiation 
related to the final state particles. 

The hard scale in the case we examine here is provided by the DM mass M > 1 TeV while the 
soft scale is the typical energy where the spectra of the final products of DM decay/annihilation 
are measured, E < 100 GeV. Even bearing in mind that weak interactions are not so weak 
at the TeV scale, one might wonder whether such "strong" electroweak effects are relevant 
for measurements with uncertainties very far from the precision reachable by ground-based 
experiments at colliders. In this context, and in view of our ignorance about the physics 
responsible for DM cross sections, it might seem that even a (9(30)% relative effect should 
have a minor impact. This is by no means the case: including electroweak corrections has a 
huge impact on the measured energy spectra from DM decay/annihilation. There are two basic 
reasons for this rather surprising result. 

• In the first place, since energy is conserved, but the total number of particles is not, 
because of electroweak radiation a small number of highly energetic particles is converted 
into a great number of low energy particles, thus enhancing the low energy ( < 100 GeV) 
part of the spectrum, which is the one of relevance from the experimental point of view. 
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• Secondly, and perhaps more importantly: since all SM particles are charged under the 
SU(2)l <8> U(1)y group, including electroweak corrections opens new channels in the final 
states which otherwise would be forbidden if such corrections are neglected. In other 
words, since electroweak corrections link all SM particles, all stable particles will be 
present in the final spectrum, independently of the primary annihilation channel consid- 
ered. 

To illustrate these facts, consider for instance a heavy DM annihilation producing an electron- 
positron pair, see Fig. 1. Clearly, as long as one does not take into account weak interactions, 
only the leptonic channel is active and no antiproton is present in the final products. However, 
at very high energies there is a probability of order unity that the positron radiates a Z or a 
W. While the spectrum of the hard positron is not much altered by virtual and real radiative 
corrections (see [17]), the Z radiation opens the hadronic channel: for instance, antiprotons are 
produced in the Z decay. Moreover, also a large number of pions are produced, which in turn 
decay to photons (tt° — > 77) and to low energy positrons (through the chain tt + — > fi + + X — > 
e + +X). At every step, energy is degraded. Because of the large multiplicity in the final states, 
the total Z energy (already smaller than the hard M scale) is distributed among a large number 
of objects, thus greatly enhancing the signal in the (10 — 100) GeV region that is measured by 
present-day experiments, like PAMELA. 

The various processes of radiation are described by fragmentation functions D EW (x, fi 2 ) that 
evolve with the energy scale /i 2 according to a set of integro-differential electroweak equations 
[18]. When a value of virtuality of the order of the weak scale \x = Mw is reached, the Z boson 
is on shell and decays. The subsequent QCD showering may be described with QCD traditional 
MonteCarlo (MC) generator tools, like PYTHIA. 

At tree level, the spectra of hard objects emerging from DM annihilation are simply propor- 
tional to 5(1 — x), where x is the fraction of center-of-mass energy carried by a given particle. 
Once electroweak corrections are switched on and 0(0,2) virtual and real corrections are cal- 
culated, the spectra -D EW (x,/i 2 = Myy) generically contains terms enhanced by log terms of 
the form a 2 In 2 M 2 /M^ and a 2 In M 2 /M^. The presence of logarithmically enhanced terms is 
well-known in the literature both in the case of electroweak interactions [14, 18] than for strong 
interactions through the Dokshitzer-Gribov-Lipatov-Altarelli-Parisi (DGLAP) equations [19] 
and is related to regions of phase space where the internal propagators become singular. In 
general the single log term is generated when two partons become collinear, while the double 
log arises when they are soft and collinear at the same time. The double log-enhanced contri- 
butions cancel in QED and QCD for physical scattering processes (Block-Nordsieck theorem 
[20]), while they are present in massive gauge theories [14]. 

The physical picture that arises is therefore the following: highly non relativistic DM par- 
ticles annihilate, producing a particle-antiparticle pair belonging to the SM spectrum. This 
pair has a very high invariant mass, therefore it stars radiating photons and gluons, but also 
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weak gauge bosons. The presence of collinear and/or infrared singularities allows to factor- 
ize leading logarithmic electroweak corrections with a probabilistic interpretation very similar 
to DGLAP equations, see Sec. 3. The exchange of virtual and emission of real electroweak 
bosons lead to the appearance in the final spectrum of all the stable SM particles, not only the 
ones initially emitted by the DM annihilation. Indeed, the higher is the mass of the DM, the 
more democratically distributed the final spectrum of DM particles is. Therefore, including 
electroweak corrections alters significantly the final spectrum of particles stemming from DM 
decay/annihilation and this has a large impact on indirect searches of DM. 

Let us close this Section by recalling that, while in this paper we only consider DM anni- 
hilation/decay to two body final states, our approach is more general and model independent. 
Indeed, the only assumptions we make are that the physics up to the DM mass scale M is 
described by the SM and that the SM may be eventually extended by interactions that pre- 
serve SU(2) L ® U(l)y gauge invariance. While these assumptions exclude cases like the ones 
considered in Refs. [10] and [21] where gauge non invariant interactions where considered 1 , a 
large number of models can be examined with the techniques we describe here. For instance, 
let us extend the SM by adding a very heavy scalar S that interacts with the SM Higgs (H) and 
leptons (L,E), through an effective operator SLEH. Then, the dominant decay of the scalar 
is a three body decay, since the two body decay S — » LE is suppressed by a relative factor 
M^y/M 2 . The framework described here applies as well, albeit with the additional complica- 
tion that the three body decay with respect to which one factorizes electroweak interactions 
provides a distribution rather than a simple 5 function. In this sense, provided assumptions 
specified above are fulfilled, our approach is completely model independent. 

3 Quantitative computation 

We now discuss in more technical terms the inclusion of EW gauge boson emission through the 
evolution equations. We start from a first principle definition of the energy spectrum for emitted 
particles and then we define the fragmentation functions as statistical objects describing the 
probability of a particle to be transformed into another with a certain momentum fraction. The 
full evolution equations for the fragmentation functions, containing EW and QCD interactions, 
are analyzed. We provide an expression that can be used to match the outcome of Monte Carlo 
codes adding EW corrections at leading order 0(0.2)- Our approach is similar in spirit to the 
one used, in a different context, in [23], with important differences. 

The relevant quantity for indirect signals of DM is the energy spectrum dNf/dx of stable SM 
particles / = {e + , e~, j,p, p, u, v,d} produced per DM decay/annihilation, where x = 2Ef/y/s 
(0 < x < 1) is the fraction of center of mass energy carried by a stable particle / with energy 

lr The analysis of the Infrared virtual corrections to gauge non invariant amplitudes, i.e. amplitudes propor- 
tionals to the higgs vev, has been recently performed in [22]. 
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Ef. For clarity we will sometimes specify the formulae assuming the case of non-relativistic 
DM annihilations, for which y/s = 2M such that x = Ef/M; it is immediate to obtain the 
corresponding formulae for DM decays, where y^i = M. 



We assume that DM initially produces two primary back-to-back SM particles, and we 
consider all relevant cases: 

^ = i e L,R e L,W ^L,R^L,Ri T L,R T L,R' ^e^e, ^ ' pP ^t^t, ^ 

qq, cc, 66, it, 77, 99, ^t,l W t,l^ ^t,l^t,l, hh} ' 

where q = u,d,s denotes a light quark; h is the Higgs boson; Left or i?ight are the possible 
fermion polarizations, and Transverse or Longitudinal are the possible polarizations of massive 
vectors, that correspond to different EW interactions. Then, the spectrum can be written as: 

^ 1 d ^™-^/+* . / = { eV-,7,P,P,,,M}, 0) 
ax cdm dm^-i ax 

with a similar formula for the case of DM decay. The "X" in this equation reminds of the 
inclusivity already discussed in Section 2. 

In each one of the possible cases /, MonteCarlo generators like PYTHIA allow to compute 
the inclusive spectrum dNf^^M, x)/dx by generating events starting from the pair / of initial 
SM particles with back-to-back momentum and energy E = M, and letting the MC to simulate 
the subsequent particle-physics evolution, taking into account decays of SM particles and their 
hadronization, as well as QCD radiation and (partially) QED radiation. 

Then, the spectra for a generic DM model that produce combinations of the two-body states 
/ can be obtained combining the various channels: 



In some DM models, primary multi-body states can be important: one can obtain the final 
spectra without running a dedicated MC code by computing the model-dependent energy spec- 
tra Dj(z) of each primary pair I (each one has energy E = zM with < z < 1) and convoluting 
them with the basic MC spectra: 

Notice that we combine particle-antiparticle pairs because we assume that they have the same 
spectra, which is true whenever the cosmological DM abundance does not carry a CP asym- 
metry. Otherwise, hadronization can be significantly affected and dedicated MC runs would 
be necessary. The indices I, J = p + p denotes a primary particle p together with its anti- 
particle p, with the same energy spectrum. Factors of two are such that for complex particles 
one has dN p /dz = dN p /dz = Dj, while for real particles (the Z, the 7, the Higgs h) one has 
dN DM ^ P /dz = 2Dj. 
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3.1 Including EW corrections 

We now come back to the basic case of DM that annihilates or decays in one primary channel / 
and discuss how to achieve the goal of this paper: obtaining a set of basic functions dNj^f/dx 
that take into account EW radiation, replacing the functions dNf^j/dx computed via Monte- 
Carlo simulations. EW radiation is a model-independent subset of the higher order corrections 
discussed above, and gives rise to specific spectra of initial SM particles, such that its effect 
can be included in the primary basic spectra by a formula similar to eq. (5): 



^(^Ejf^M^^, (6, 

where Df^j(z) is the EW I — > J EW parton distribution: the J spectrum produced by initial 
I. Our normalization is such that we have the uniform normalization Df^j(z) = 5jj5(1 — z) 
at tree level for both real and complex particles. Some comments are in order: 

i) When including higher order effects, one must avoid overcounting and take into account 
that MC codes already include some particularly relevant higher-order effects: showering pro- 
duced by strong (QCD) and electromagnetic (QED) interactions, up to details. 2 

ii) For initial particles that do have strong interactions, eq. (6) misses the interplay between 
EW and QCD radiation; this limitation is not a problem because, as expected, in such cases 
EW radiation will turn out to be subdominant with respect to QCD radiation. 

iii) For initial particles that do not have strong interactions, eq. (6) holds at leading order in 
the weak couplings: first they must do an EW splitting, and next one can add QCD splittings 
neglecting EW radiation. We emphasize an important different between e.g. a Z — > qq splitting 
and the same Z — > qq decay: in the decay the invariant mass of the qq pair is equal to the Z 
mass (such that Z — > ti is forbidden by the heaviness of the top t); in the splitting the invariant 
mass can be much higher, and we approximate it as zM. This higher invariant mass strongly 
affects the subsequent QCD radiation from quarks, which is more abundant in the splitting 
case, leading to a higher multiplicity of p and 7. 

In Appendix A we give a detailed discussion of the interplay between EW and QCD radiation 
and of the level of approximation introduced by using eq. (6). 



x 



3.2 Computing the EW parton distributions 

We define Df^j(z, fi 2 ) as the probability for a given parent particle I with virtuality of the 
order of /1 to become a particle J with a fraction z of the parent particle's energy mediated by 

2 We must include via eq. (6) only all those effects not included in MC codes. Existing MC codes have their 
own peculiarities, e.g. Phytia automatically includes 7 radiation from charged particles but not from the W . 
Of course, an alternative more precise approach, that we do not purse, would be implementing the missing EW 
radiation effects into some existing MC code. 
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EW interactions. At large virtuality, they take the tree level values: 

DfZ(z,» 2 = s) = SjjSil-z); (7) 

At low virtuality /i 2 ~ M^, they are the functions we need: Df*j(z) = Dj™j(z,fi 2 = M^). 
The evolution in the virtuality is described by integro-differential equations, that involve a set 
of kernels 3 Ppj^j(x, fi 2 ) that have been derived in [18]: 

dD ^ z / ] =-SE / v p^kM d^Az/v^ 2 )- (8) 

k 

Since we work at leading order in the EW couplings, eq. (8) with the boundary conditions of 
eq. (7) is solved by: 



a 2 f s d/j 2 
2tt Art /i 2 



Differently from QED and QCD, the EW kernels feature infrared singular terms proportional 
to ln/i 2 , so that the solutions (9) also include double logs beside the customary single logs of 
collinear origin: 

DfZ(z) = D 2 (z) In 2 + D x {z) In + D (z). (10) 

IV1\Y lV±w 

Our goal is to include the model-independent logarithmically enhanced terms. Electroweak 
radiation from the initial DM state is of course model-dependent: since DM is non-relativistic 
this effect only contributes to the non-enhanced terms Do, that we neglect. Notice however 
that for our purposes we need to include terms of the form (lnx)/x that are relevant in the 
region x — > 0; this is discussed in detail in subsection 3.4. 

3.3 Splitting functions 

The leading-order parton distributions Dj^j(z) can be computed by using the partonic splitting 
functions P summing over all possible SM splittings [18]; the relevant splitting functions are 
here collected in table 1. 

A concrete simple example allows to clarify the procedure and the normalization factors: we 
consider DM producing an initial generic Fermion-antiFermion pair with FF invariant mass 
a/s ^> mp. We assume that F has charge qp under a generic vector V with mass My <C 
and gauge coupling a. (In the SM the vector could e.g. be a Z and the fermion a neutrino). 
The splitting process F —> FV gives rise to: 



D M {z) = 6(1 - z) 



1 + ^P^ F 



+ ^W04 D F ^y{z) = %^P F ^ V {z). (11) 

Z7T Z7T 



3 In this work we indicate with P(x, fj, 2 ) the unintegrated kernels, while the splitting functions P{x), obtained 
by integrating in /i 2 , depend only on the energy fraction x; a list of the relevant splitting functions is given in 
Table 1. 
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splitting 1 — > x + x' 



splitting function: real and virtual 



Fq,m — > Fq,m + Vm Pf^f 
Fqm Vm + F M Pf^v 



Vm y Vm + Vm P 



Vm ^ Vm + Vq P'y^y 



Vm ^ K, + V M 

F — »■ F + S 

F — >■ S + F 

5^ F + F 



Fy-> 7 
pYuk 



P 



F^-F 
Yuk 



P 



Yuk 
S->F 



1 + X z 



1 — x 
l + (l-x) 2 



pvir 
r F^F 



V 


F + F 


Py^F — 


[x 2 + (l-x] 


2 Y 


5"m — >■ 


5"m + Vm 


Ps->s — 


2 X L(l 


- x) 








1 — X 




Sm — >■ 


Vm + Sm 








Ps^v = 


2 1 ~ X L(x) 










X 




v -> 


S + S' 


Py->s — 


x(l — x)£ 





L{x) Ppfj^y 



pvir 
r V-+F 



X 


1 




X 




X 




1 




X 




X 




1 




X 



2 
2 
2 

(i-x)e 

x£ 



L{l-x) 



pvir 



pvir 



pvir 

1 — X 

X 



U _ I 1 

3£ _ e 

~2 ~ 2" 

2£ 



2£- 
2£ 



6 

L(s) + x(l - x)£ 



X 



L(x) +x(l-x)£ 



L(l-x) + - — -£ + x(l - x)£ 



x 

3>Yuk,vir 
F^F 

3>Yuk,vir 

3>Yuk,vir 
S-+F 



2 

£ 

2 



Table 1: Generalized splitting functions for massive partons. V denotes a vector, F a fermion 
and S a scalar; Vm denotes a vector with mass My, Vq a massless vector, etc. The function 
L(x) is defined in eq. (18) and £ = Ins/My. 



By replacing F — > S one obtains the corresponding result for a pair for Scalars, and so on. 

The first term describes virtual corrections arising from one-loop diagrams, and the second 
term describes real emission corrections. We define: 



pvir 



dz P, 



for any I and J; e.g. from Py^y in table 1 we have: 

3 



pvir 
r V-*V 



where 



In 



My 



(12) 



(13) 



While both virtual corrections, related to P/_">j, and real corrections, related to P/_*.j, are listed 
in table 1, the simple relation (12) holds between them. The relationship between real and 
virtual contributions is dictated by the unitarity of the theory (see e.g. [17] for a more detailed 
discussion). Intuitively, this just amounts to say that when an F radiates, it disappears from 
the initial state. 
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One can verify that the partonic distributions Di^j satisfy a set of conservations laws (with 
corresponding identities for the splitting functions P) 



The conservation of splitting probability: 



2 
j 



dxD\% + / dxD^j = 0, (14) 
o Jo 



where the factor 1/2 accounts for the fact that one particle splits in two. 

• The conservation of total momentum: 

V I dx x Di^j = 1. (15) 
j Jo 

• The conservation of electrical charge: 

V l dxQj Dj^j = Qj. (16) 
j Jo 

By combining these splitting functions with the appropriate electroweak couplings (including 
the top quark Yukawa coupling) we get the electroweak splittings among SM particles, described 
by the Dj^j functions, explicitly listed in Appendix C. It is a simple exercise to verify that 
conservation laws are satisfied. For completeness, we also list the splittings involving photons 
(to be dropped if already included by MC codes), assuming for simplicity a photon with mass 
M w ; sending it to zero gives infrared divergences that can be regulated and dealt with using 
well known techniques that we do not need to discuss here. 

3.4 Splitting functions for massive partons 

We see from table 1 that PJ^j is explicitly given by linear or quadratic polynomials in In s/My\ 
the vector mass provides an infra-red regulator. This is unlike in the standard application of 
partonic techniques, where all partons are massless and some infra-red regularization is needed. 

In this subsection we briefly describe how we generalized partonic functions to massive 
partons (such as the W, Z, t in the SM), and why particle masses modify the kinematics affecting 
the log-enhanced terms, as encoded in a non trivial universal function L. 

Let us consider the splitting process i — >■/ + /' of a particle i into massive partons /, /': 
the corresponding splitting functions can be non zero only in the kinematically allowed ranges: 

rrif mt' , „ 

-4- < * < 1 - -f , (17) 

where x = Ef/Ei, and similarly for x' = 1 — x = Ef/Ei. The modification due to nonzero 
masses in the kinematical range is power suppressed and thereby mostly negligible; however 
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Z from W radiation, M = 3 TeV 




full 3 body 



0.1 



0.01 



0.03 



0.1 



0.3 



x = E/M 



Figure 2: Splitting function between massive vectors. Blue curve: naive result. Dashed 
curve: our result. Red curve: full 3-body result in the Minimal Dark Matter model. 

the small forbidden regions at small x, x' are relevant when vectors are emitted, because the 
corresponding splitting functions have 1/x and 1/x' = 1/(1 — x) divergences, see table 1. 
Thereby small parton masses must strongly affect such splitting functions, replacing the usual 
1/x and 1/x 1 divergent terms with new functions that vanish at the kinematical boundary. 

We therefore need to correct the expressions for the splitting functions that are singular in 
the limits x — > and x' — > (soft / and soft /'). These singularities arises when soft vectors 
are emitted, such that the correct expressions can be obtained using the eikonal approximation, 
which dictates the amplitude for emitting a soft vector. Since we stop at leading order, we can 
precisely define the splitting functions in terms of the energy spectra of the particles produced 
in three body scatterings: we just need to integrate such squared amplitude over the massive 
phase space in the region that produces the leading singularities. Broadly speaking, what 
happens is that the upper and lower limits of integration that appear in (9) are modified and 
become x-dependent; computations and checks are performed in Appendix B. The resulting 
splitting functions for massive partons are listed in table 1, and, for x — Y 0, they do not depend 
on spin, e.g. Ps^v — Pf^v — Pv^v, as dictated by the eikonal amplitude. They contain the 
universal kinematical function L(x), that replaces the usual t = las /My that holds for massive 
partons: 



This function indeed vanishes below the kinematical threshold (x < 2Mv/y/~s) and reduces to 
In sx 2 /M v well above it. This means that small parton masses, apart from providing kinematical 
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2.06 
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DM DM -> W^W^t 


15.1 


18.9 


24.8 


31.3 


39.8 


53.2 


1.50 


1.92 


2.62 


DM DM -> W£W£ 


15.5 


19.6 


25.8 


32.3 


41.6 


55.9 


1.52 


1.94 


2.66 


DM DM /i/i 


27.3 


30.5 


35.5 


59.9 


67.1 


78.6 


2.22 


2.60 


3.26 


DM DM gg 


22.6 


34.6 


48.9 


47.5 


73.9 


107. 


2.47 


3.89 


5.71 



Table 2: Total number of e + , 7, p for a few DM annihilation channels. 



thresholds, give rise to extra \nx factors with respect to the standard case of massless partons, 
which become numerically relevant at small x, x'. 

This was not noticed before, and Fig. 2 exemplifies its relevance. The dashed curve in Fig. 
2 shows the splitting function between massive vectors (e.g. relevant for Z radiation from W ± ): 
it significantly differs from the massless splitting function (upper curve) even away from the 
kinematical boundaries, and it closely agrees with the full result of a 3-body computation in a 
specific model, which also includes not log-enhanced terms. These terms are subleading in the 
present example, where we here considered DM DM — > W^W^f annihilations with M = 3 TeV. 

4 Results 

Our main results are the energy spectra of all stable final SM particles / from any two-body 
DM non-relativistic annihilations or decays. We will make the code freely available in [24], and 
show here some examples. 

In table 2 we show the total number of e + , p and 7, including EW radiation, for a few values 
of the DM mass M. Without EW radiation, increasing M would just boost all particles, such 
that dN/ dx and consequently N does not depend on M. With EW radiation, the total number 
of particles increases with the DM mass: the typical effect is about a factor of 2 going from 0.3 
to 3 TeV. Of course, such dependence on M is similar to the one already present in the quark 
channels (bottom row), due to QCD radiation. 

Coming to the energy spectra, these new particles appear at low energy. In Fig. 3 we 
consider annihilations of DM with mass M = 3 TeV and compare the spectra with (continuous 
curves) and without (dashed curves) EW corrections. 

• In the top row we consider DM annihilations into W + W~ with transverse (left) or longitu- 
dinal (right) polarization. The final spectra do not significantly depend on the polarization 
of the W, but this is an accident: a Transverse W splits into all quarks as demanded 
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Figure 3: Comparison between spectra with (continuous lines) and without EW corrections 
(dashed). We show the following final states: e + (green), p (blue), 7 (red), v = [y e + + z/ T )/3 
(black). 
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Figure 4: DM signals in the e + (left) and p (right) fraction, with (dashed) and without (dot- 
dashed) electroweak corrections for two DM models that can fit the PAMELA e + excess: Minimal 
Dark Matter (upper) or a muonic channel (lower). The gray area is the predicted astrophysical 
background and the red area is the prediction adding the full DM contribution. 



by their weak interactions, while a Longitudinal W (being the charged Goldstone in the 
Higgs doublet) mainly splits into t and b quarks, as demanded by the top quark Yukawa 
interaction. In both cases EW splittings increase the total number of final particles (e + , p, 
7, v) increases by a factor of almost 2. The biggest effect is in the p spectra: without EW 
corrections they are strongly suppressed below E p < m p ■ (M/Mw) ~ 100 GeV, because 
of the boost factor M/Mw of the W. Adding EW corrections, a W splits into quarks 
with lower energy, producing p also at lower energy. This is important for interpreta- 
tions of the PAMELA e + excess in terms of annihilations of very heavy DM particles, 
e.g. M = lOTeV as predicted by Minimal Dark Matter [25]. EW corrections make more 
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difficult to avoid unseen effects in p by having all the DM effects above 100 GeV where 
we do not have data yet. As illustrated in Fig. 4, in view of EW corrections there is some 
tension with data, even assuming the MIN model of diffusion of charged particles in the 
galaxy [26] (which gives the minimal amount of p compatibly with cosmic rays data and 
theory) . 

• We do not plot our result for DM annihilations into the Higgs boson, because its mass 
and decay modes are not yet known. The Higgs channel is similarly affected by EW 
corrections as the Wl and Zl channels (because they are the Goldstone components of 
the Higgs doublet), and the main effect is again the h — > t splitting induced by the top 
Yukawa coupling. 

• In the middle row of Fig. 3 we consider DM annihilations into charged leptons: e^e^ 
(left) and (right). The spectra are significantly affected by EW corrections, because 
leptons split into W, Z bosons, finally producing quarks, and consequently a copious tail 
of e + ,7, p at lower energy. The main new qualitative feature is the appearance of p 
from leptons, and in Fig. 4 (lower row) we show that DM annihilations into with 
M = 2 TeV (a scenario motivated by the PAMELA e + and FERMI e + + e~ anomalies) 
also gives a possibly detectable excess of p. We here assumed the favored MED model 
of propagation of charged particles in the Milky Way [26]; p can be suppressed down to 
a negligible by considering the MIN model and/or DM annihilations into /U^u^- Indeed 
EW effects are more significant for left-handed leptons which have SU(2)l interactions 
than for right-handed leptons which only have U(l)y interactions. 

• The bottom left panel of Fig. 3 shows that, in view of EW interactions, DM annihilations 
into neutrinos also induce a significant spectrum of 7 and some e + , p. 

• The bottom right panel of Fig. 3 shows how EW corrections affect DM annihilations into 
77: an hypothetical 7 line at high energy E 1 = M ~ few TeV must be accompanied by a 
comparable flux of 7 with lower energy ~ (10 — 100) GeV, where we have more data. 4 

• Finally, DM annihilations and decays into quarks are negligibly affected, so that we do 
not show our results. 

The same results hold for DM decays. EW corrections are irrelevant in models where DM 
decays or annihilates into hypothetical 'dark' particles, lighter than the weak scale, that finally 
decay back to SM particles. 

There is a striking difference between the two examples of e + fraction in Fig. 4: EW 
corrections are relevant in the upper case (Wt channel) and negligible in the lower case (//^ 
channel). At first sight, this is surprising, because Fig. 3 shows that in both cases EW 

4 This effect, already partly present in QED, seems to be not implemented into MonteCarlo codes. 
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corrections induce significant low energy tails to the e + spectra at production. The qualitative 
difference is due to e + energy losses in the Milky Way that also gives a tail of e + at low energy: 
we explain how to understand which effect is dominant. The e + flux can be approximatively 
computed neglecting galactic diffusion and taking into account only energy losses, with the 
result 

dn e - _ dn e + 3m 2 e 1 ov fP®\ 2 f M dE , dNe (lq s 
dE ~ dE ~ Aa TU(B E 2 2 \MJ J E dE' [ ' 

where <jt is the Thompson cross section, u e is the energy density in radiation and magnetic 
fields and p e is the DM density, both at the location of the solar system. This means that 
the e + flux at a given energy E < M is proportional to the number of e + produced with 
energy between E and M. In the case of DM annihilations into leptonic channels, the tail of 
e at lower energy produced by EW radiation contains a number or at most comparable to 
the amount of e ± already present at E ~ M: therefore EW corrections negligibly affect the 
positron fraction. This is not the case for DM annihilations into W ± , where instead the tail at 
low energy is the dominant component of the total e number, such that EW corrections more 
significantly affect the e + fraction. 



5 Conclusions 

In this paper we computed the energy spectrum of SM particles stemming out from DM annihi- 
lation/decay. We have shown that EW corrections have a relevant impact on such spectra when 
the mass M of the DM particles is larger than the EW scale Mw- Soft EW boson emission is 
enhanced in the collinear and infrared regime and this leads to In 2 M 2 /M^ enhancement fac- 
tors. The result of the inclusion of EW corrections is that all stable particles are present in the 
final spectrum, independently of the primary annihilation/decay channel. For instance, even if 
the DM particles annihilate/decay into light neutrinos, EW corrections cause the appearance 
of hadrons and photons in the very final spectrum. The inclusion of EW corrections is therefore 
an essential ingredient in order to have a physical picture of the correlated energy spectra of 
final stable particles. Our quantitative results may be inferred from Fig. 3 where the energy 
spectra dN/dE of stable SM particles are presented with and without EW corrections. These 
spectra are the necessary ingredients to predict the flux for indirect searches once the effect 
of diffusion and galactic energy loss are included. As a rule of thumb we may say that EW 
corrections are important in determining the final flux of stable particles whenever 

• the final flux of stable particles is dominated by the low energy tail of the dN/dE. One 
example is the case of DM annihilation/decay into gauge bosons and e final states (see 
fig. 3). This point becomes more relevant in the present experimental situation, where 
we mostly observe particles below 100 GeV, possibly much below the DM mass. 
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• the final flux of stable particles is absent when EW corrections are not taken into account. 
One example is the case of DM annihilation/decay into leptons and antiprotons p final 
states (see fig. 4). This point is important also for neutrino fluxes from DM annihilations 
in the sun or in the earth, because all SM particles, even those that loose energy in 
matter before decaying into neutrinos, can radiate a W or a Z that promptly decays into 
neutrinos. 

EW corrections may also significantly affect the fluxes of particles generated from heavy grav- 
itino decays in supersymmetric theories. The computation of these energy fluxes is crucial in 
studying the dissociation of the light element abundances generated during a period primordial 
nucleosynthesis. 

We expect that EW radiative corrections have a minor effect on the freeze-out cosmological 
DM abundance, because it is determined by just the total non-relativist ic annihilation DM cross 
section. In the case of energy spectra instead, as explained in this paper, the low energy tails 
can be enhanced by orders of magnitude, while the high energy part of the spectrum is mildly 
depleted. The net effect on the total number of final particles typically is an enhancement by a 
factor of 2. We conclude that, when DM is around or above the TeV scale, one must take into 
account radiative EW corrections. 

We computed EW corrections at leading order. Although we cannot give a sharp answer, 
we argued that it is not necessary to resum higher order EW corrections as long as DM is 
not too heavy: CK2 ln 2 (M 2 /M^)/27r < 1 at M < 100 TeV. Indeed, higher order corrections 
are expected to give small effects for total cross sections at the TeV scale: a one loop effect 
of the order of 30% means that one expects higher order effects to be at the 1% level. It is 
more difficult to guess how low-energy tails might be affected without performing a dedicated 
computation, which could be done implementing EW corrections in a MonteCarlo: we provided 
splitting functions for massive partons and all ingredients. 

Acknowledgments: D.C. during this work was partially supported by the EU FP6 Marie Curie 
Research & Training Network " UniverseNet" (MRTN-CT-2006-035863). P.C. wishes to thank Luigi 
Lanzolla for many useful contacts and discussions, without which his work would have been much 
more difficult. We also thank G.F. Giudice for useful conversations. We used the high-statistics 
Pythia spectra computed on the Baltic Grid by Mario Kadastik, to appear in [24]; this work was 
supported by the ESF Grant 8499 and by the MTT8 project. 

A Evolution Equations 

At present the energy spectrum dNf^^/dx are computed with MC generators like PYTHIA by gener- 
ating events starting from the pair I of initial SM particles with back-to-back momentum and energy 
E = M, and letting the MC to simulate the subsequent particle-physics evolution, taking into account 
decays of SM particles and their hadronization, as well as QED and QCD radiation. As it is evident 
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the MC output are nothing else that the full QCD+QED fragmentation functions D^_^j (x) that 
can be identified with dNj^j/dx. The fragmentation functions are related to a nonperturbative as- 
pect of QCD, so that they cannot be precisely calculated by theoretical methods at this stage. The 
situation is similar to the determination of the PDFs, where high-energy experimental data are used 
for their determination instead of theoretical calculations. The fi 2 evolution for the fragmentation 
functions is calculated in the same way as the one for the PDFs by using the timelike DGLAP evo- 
lution equations. The splitting functions are the same in the Leading Order (LO) evolution of the 
PDFs; however, they are different in the Next-to-Leading Order (NLO). Explicit forms of the splitting 
functions are provided in [18]. The evolution equations are then essentially the same as the PDF case. 

The correct determination of the energy spectrum dNj^j/dx of the final stable particle /, ob- 
tained through the particle-physical evolution of the initial pair / of SM particles with back-to-back 
momentum and energy E = M, needs the solution of a full set of evolution equations, including both 
strong and electroweak interactions. 

In previous works the QCD DGLAP formalism has been extended to EW interactions [18]. The anal- 
ysis of mass singularities in a spontaneously broken gauge theory like the electroweak sector of the 
Standard Model has many interesting features. To begin with, initial states like electrons and protons 
carry nonabelian (isospin) charges; this feature causes the very existence of double logs, i.e. the lack 
of cancellations of virtual corrections with real emission in inclusive observables [14]. Secondly, initial 
states that are mass eigenstates are not necessarily gauge eigenstates; this causes some interesting 
mixing phenomena analyzed in [16]. 

On a general setting, the mathematical structure of the full system of EW and QCD evolution 
is provided by the following set of integral differential equations with kernels pQ CD / EW f or different 
gauge boson exchange (we omit for the moment the explicit QED evolution equations included in the 
EW part of the SM and the various index parametrizing the flavour and quantum numbers): 5 

fi 2 ^D(x,/j 2 ) = — D®P® CD d(A< fi< M) + °^ D®P EW 6{M W <^<M) 

= D® (^ P QCD + ^pEw\ e{M w<l x<M) 
\ In 2tt J 

+ ^ D®P® CD 9(A Q cd<V<M w ), (A.l) 

where we have made explicit the running range: for the EW interactions from M to M\y where the 
gauge boson mass is freezing the running and for QCD corrections from M to Aqcd where the non- 
perturbative effects generate an effective cutoff to gluon exchanges. In practice, from Aqcd to Myy 
the running is purely dictated by QCD (remember that QED interactions are not shown for simplicity, 
in this case the running scale for photons is until the electron mass m e ) while the EW-QCD interplay 
starts only above the M\y scale. Note also that forgetting the underlying index structure can bring 
to the wrong conclusion that, simply due to the fact that a s > at2, the EW part is just a correction to 
the QCD dominant part. This would be a mistake as in many interesting annihilation channels (like 
all the leptonic ones or the W^, Z or h) the QCD part is simple zero. 

5 The precise index (flavour) structure is given in the Appendix C, the generic index structure is of the 
form y^-J^Di^j = ^ J2 K Di^K ® Pr^j while the <g)-operator means (/ <8> g){x) = f(z) ® g(x/z) = 
f x dz/zf(z)g (x/z) = f dz ^ dyf(y)g(z)S{x - zy). 
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A numerical solution to the full (EW+QCD) problem is of course out of reach. Nevertheless, 
we can find some reasonable approximation taking advantage of the fact that we can simulate the 
pure QCD evolution also in the non perturbative regime with MC codes and that EW theory is in 
the perturbative regime. Technically specking, the evolution equations are Schrodinger-like equations 
with a time dependent Hamiltonian, where time is replaced by the /i 2 variable and the Hamiltonian 
by the P- kernel. A formal solution can be parametrized with the evolution operator: 



u(z,i4,i4) 




(A.2) 



where V^i is the /^-ordering operator and I the identity in the flavour space. Due to the linearity of 
the Eqs. (A.l) 6 we can then formally write the full solution as: 



D(x,M 2 ,A 2 QCD ) = U(z,M 2 ,M 2 v ) 



Wi ll QCD 



(A.3) 



where we have separated the running from M to Myy inside the evolution operator U (that can be 
perturbatively expanded as soon as a s In M 2 /M^ and a 2 In 2 M 2 /M$ v are smaller than unity) and the 
purely QCD piece D^ CD (x, M^, Aq CD ) encoding also the non-perturbative low energy physics. In 
order to keep under control further simplifications we need also to know the matrix flavour structure. 
We display it under the form of a simplified four-dimensional space spanned by Z=leptons, W = 
(W ± ,Z,j), g=quarks and g=gluons, for the EW and QCD kernels, pQ CD and P EW : 
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(A.4) 



The above matrices do not commute; the EW and QCD sectors are connected through the channels 
W — > q, q —7- W and q —> q, furthermore the leptonic and hadronic sectors are connected through 
the mixed W — > q, I, I, q — )• W channels. Reasonable approximate solutions are related both to 
the possibility to expand perturbatively the general solution (B.15), and to the outcome from MC 
generators which take automatically into account the full QCD plus QED evolution from M to Aqcd 
scales (m e for QED). 

One way of proceeding is to define pure EW (£> EW ) and QCD (D QCD ) fragmentation functions 
that evolve with their respective kernels, see eq. (A.4), for the energy range Myy < n 2 < M 2 : 



, 2 ^- 2 D^( X ,, 2 ) 
whose formal solutions are: 



^£ £)QCD 
2vr 



PQ° D and 



2vr 



P 



EW 



D EW ^ CD (x,M 2 ,M^) 



P M 2 e 



rM 2 dp? 

Jm w 



pEW/QCD^l 



I. 



(A.5) 



(A.6) 



3 It might be useful to remember the property D(z, M 2 , Myy) 



D{x/z,M^,A 2 QCn ) 



D(x,M 2 ,A 



QCD) 
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Then we introduce a new factorized EW ® QCD fragmentation function: 



D{x^ 2 ) = {D 



EW 



^QCD )( 



with 9(M W < fx<M). 



(A.7) 



This is clearly not a solution of the true evolution equations (A.l) but can be a useful approximate 
solution. In order to relate the true solution D of eq. (A. 3) with the new function D satisfying eq. 
(A.7)), we can use the fact that in the M w < fj, 2 < M 2 interval we are in perturbative regime also for 
the QCD side. Knowing that for two generic non-commuting operators A and B: 

1, 
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;[A,B] + 



e A e B 



(A. 



and identifying A = J } 



M 2 dfi 2 (a 2 
M 2 — ~~ 
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M 2 djj? (a 



M 



W ^ 



i^n -P^ CD ), we can approximate the evo- 



lution operator U at any order in a s ,2- In particular at second order in 0(a s2 , 0^02) we have: 

a s ot2 
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D EW ^ D QCD 
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M 2 W F 
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(x,M 2 ,M w ), 



where: 



[pQ CD , p E % = (pQ° D ® p ew - p ew ® pQ° D ) 



(A.9) 



(A.10) 



and the • • • stand for the fact that there is an infinite series of commutators with coefficients of order 
a 2 n+1 with m + n > 1. Starting from this expression we can write the perturbative relation 

between the exact solution D and the present outcome of the MC codes Z)Q CD : 
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(A.ll) 



The first order corrections to such a formula are obtained expanding also D EW at one loop: 

dfi' 2 
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P EW (x,v 



12 \ 



/<-' 



where we have explicitly shown the arguments x and \i of the matrix p EW . 

At order 0(a<i In 2 M 2 /M w ,a s In M 2 /M w ) we find that the energy spectrum for the process / 

can be therefore written as in eq. (6): 



(A.12) 
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V xi 



dx 



(A.13) 



where in the last passage, to be consistent, we have written EW in order to stress that only massive 
and Z are included while QED is already encoded in the MC. 
Through this expression one can match the MC code with the first order EW corrections. 
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\ W + ( Pl ) W + ( Pl ) W + ( Pl )/ 

Figure 5: Soft gauge boson real emission from spin-1 particle. It can be read as the sum of 
three scalar currents. Considering the process y/s — > p\ + k + P2 we show explicitly only the 
bremsstrahlung contribution from the p\ final leg. 



B Eikonal approximation and the improved splitting func- 
tions 

The standard partonic approximation holds in QED and in QCD for the emission of soft massless gauge 
bosons (photons or gluons) from partons, showing the presence of universal logarithmical factors of 
collinear origin that multiply the usual splitting functions. This approach can't be naively applied to 
the electroweak case when a massive gauge boson, such as the W, is involved in the splitting process; 
considering for defmiteness i — > f + /' and defying x = Ef/Ei, in fact, the allowed kinematical range 
for this latter is: 

mt mf ,„ . 

~^<x<l--f, (B.l) 

where particle masses act as cut-off for the soft singularities at x — ¥ 0, 1. These boundary regions in 
(B.l) are therefore extremely important and the standard partonic approximation have to be improved, 
introducing extra lnx and ln(l — x) terms, well justified by the kinematical proprieties of the splitting 
process. 

In this Appendix we derive our improved splitting functions for massive partons, following the logic 
outlined in section 3. In B.l we use the eikonal approximation, that describes the amplitudes with soft 
gauge boson emission. We integrate it over the phase space using for it the Sudakov approximation in 
B.2 and using the exact expression in B.3: the Sudakov parametrization, commonly used in literature, 
don't respect the boundaries in (B.l). In B.4 we introduce, through an explicit example, the collinear 
approximation and its propriety of factorization. Finally, in B.5 we compare our results with those of 
a full full three body calculation (exact amplitude integrated over the exact phase space). 

B.l The eikonal amplitude 

As well known, the spin of the emitting particles (scalar, fermion and vector) becomes irrelevant in 
the eikonal limit: for defmiteness, and without losing generality we here consider the real emission of 
a particle with momentum k described through the three gauge boson vertex 3g. 
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Considering Fig. 5 and using the conservation of momenta in the splitting vertex p = p\ + k we 
can write: 

iM " { 2^k \- £ * {k) ' £ * {Pl){k ~ Pl)p + £ > l)2pi " £ * {k) ~ £ *P {k)2k " £ * (pi) ] } ' (B - 2) 

where we have indicated with iAi p the remaining part of the amplitude and taken for simplicity k 2 = 0. 
Roughly speaking using the eikonal approximation it is possible to neglect the soft momenta in the 
numerator in front of the hard one and in this example we can discuss the two opposite situation in 
which either p\ or k are soft. 

Considering Ward identities the first term in (B.2) vanishes in the case of transverse gauge bosons, 
while for longitudinal degrees of freedom we refer the interested reader to [27] for an analysis of the 
electroweak symmetry breaking effects; therefore for our purpose we have: 



where: 



iMP 1 ~p^k \- £ > l)pi " £ * {k) ~ £ P {k)k " £ * (pi) ] f ' (B - 3) 



the first term reconstructs the hard scattering amplitude iAi p £*(px) and it survives when k is 
soft; 

the second term reconstructs the hard scattering amplitude iM p e* p {k) and it survives when p\ 
is soft. 



Squaring the amplitude (B.3), we sum over polarizations using the axial gauge 

p 2 ■ k 

and as a result the eikonal limit leads to the factorization of the process in the product of the hard 
cross section times an emission factor integrated over the allowed phase space of the soft particle: 



= -r + k " P2 + ] f p \ (b.4) 



T( \ 2 2k -P2 d 4 p! 2 | 

Z{PI) = 9 ( Pl ■ k)( P2 ■ Pl ) (2^ 6 ^M>o ^ lSoft ' 

m = 9 \p 1 %llk ) Wf 5{k2) ^ ksofL (R6) 

The two integrals can be obtained through the exchange pi f> fc. 

B.2 The Sudakov parametrization 

We consider now the explicit evaluation of the eikonal integral in (B.5), and in order to perform this 
calculation we choose a convenient parametrization of the external momenta; fixing two basis vectors: 

P = (E,0,0,E), P = (E,0,0,-E), with: s = 4E 2 ~ AM 2 , (B.7) 
the Sudakov parametrization consists in the following decomposition for the soft momentum p\ : 

Pl = X P + xP - k± = (E(x + x),-k t ,0,E(x -x)) , (B.8) 
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where, without losing generality, we have taken the spatial component of k± along the x direction. 

In order to highlight in a simple way the logarithmical behavior of the eikonal integral we choose 
to work considering as first step the massless case and approximating the two hard momenta k — P 
and P2 — P, such that 



d A pi = S -^dk^dxdx. (B.9) 



The eikonal integral takes the form: 



X(pi) = dk? dxdxS (sxx — k?) — — , (B.10) 

IT XX 

and logarithmical singularities clearly arise in two opposite kinematical regions: 

• x 3> x: the soft gauge boson p\ is emitted along the k direction; integrating over x using 
x = kf /sx, the condition i>i becomes an upper bound for the transverse momentum k 2 <C sx 2 
and therefore in terms of fragmentation functions in the p\ soft limit x — > we obtain: 

A^O = r dkl 02 2 1 = 02 2 In (B.ll) 

that vanishes when x = Mw / Vs. 

• x -C x: integrating over x using the relation x = k 2 /sx, the eikonal integral gives exactly the 
same logarithmical result previously discussed but in an opposite kinematical configuration since 
the soft gauge boson p\ is emitted now along p2 direction. 

In order to clarify the consequences of the symmetry p\ -R- k between the two eikonal integrals, we 
can now discuss in more details the Sudakov parametrization in the region x>i, First we generalize 
eq. (B.8) writing for k and p2- 

k = zP + zP + k ± = {E{z + z),k t ,0,E{z-z)), P2 = yP = {yE,0,0,-yE), (B.12) 

and using the on-shell conditions in order to eliminate z, x writing z = k 2 /sz, x = k 2 /sx. The 
conservation of energy and spatial momentum gives the following relations between the kinematical 
variables x, y, z: 

and we can therefore generalize the result for k soft just considering the substitution p\ — > k ==^ 
x — > z = 1 — x, and as a consequence the kinematical end-points for the x variable in the Sudakov 
parametrization are: 

^<*<1-^. (B.14) 



Comparing eq. (B.ll) with: 
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where P{x, k 2 ) is the usual unintegrated splitting function, we obtain its leading behavior in corre- 
spondence of the two kinematical limit in the Sudakov parametrization: 

^^f: ^Sud ~ - x L{x)\ SnA , (B.16) 

x-H-^J: ^Sud ~ Y^- x L{l-x)\ Sud , (B.17) 



with: 



2 

^)lsud = lnv7 2 -. (B.18) 

1V1 W 



B.2.1 Parton masses and the lower limit of integration 

The upper bound on the integration over k 2 is dictated by the kinematical proprieties of the collinear 
emission, and can be studied even working in the massless case; parton masses, contrarily, affect in 
a relevant way the lower bound of integration. To discuss this point, we need first to generalize the 
Sudakov parametrization in eq. (B.12): it's straightforward to verify that one can take into account 
of the on-shell conditions k 2 = m 2 and p\ = m\ just redefining k 2 — >■ k 2 + rrvl for k and k 2 — > k 2 + m\ 
for p±. Then from the propagator of the collinear emission p — > p\ + k we have: 

1 x(l — x) 



(pi + k) 2 — m 2 k 2 +m\ + x{jn\ — m 2 ) — m 2 x(l — x) 
Depending on which particles are massive, three different situations arise. 



(B.19) 



1. Only the emitted vector is massive. This happens e.g. in EW interactions of the W, Z vectors. 
Assuming vn\ = M^, m 2 = m 2 = 0, eq. (B.19) becomes: 

1 _ X(l-X) X(l-X) 2 2 

Jj^+W - k 2 + xM 2 w kj~ m XMwh (B - 20) 

where the latter passage holds in logarithmical accuracy. Integrating over k 2 we have for the 
x —7- 1 singularity of the EW splitting function Pp-+p: 

ft^.-l-taja^! _,),„, (B.21) 
and therefore the lower x-dependence don't affect the soft limit x — > 1. 



2. Two particles are massive. This happens e.g. in the electromagnetic coupling of the W: m, 



2 



in 



I = Myy, m i = 0; eq. (B.19) becomes: 



{pi + k 

Integrating over k 2 we have for the x — > 0, 1 singularities of the splitting function Py^\ 



2 , s 1 -x) 2 2 , sx 2 2 r , Nl 2, s 
Py^-Y In — ^ + - In — - ^(1 - z)L,h + - In — (B.23) 

^ 7 1-x x^ 2 , x x 2 M^ 1-x V 7ISud x M 2 / V 7 

and the lower x-dependence affects the x — > singularity for the soft photon. 
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3. All three particles are massive. This happens in the massive three gauge boson vertex, m 2 = 
m\ = m\ = Myf \ eq. (B.19) becomes: 

, rrA M 2 « - (1 - x + x 2 )M^]. (B.24) 

(pi + k) 2 - Mfo kf 

Integrating over kf we have for the x — > 0, 1 singularities of the splitting function Py^y. 

2 s (l — x) 2 2 sx 2 2 2 
P V^V~ z -In- +- ln 7, ,, ■ „2Mi,r2 -i Z i ( 1 - a 0lsud + - L ( a 0lsud» 



1 - x- (l-x + x 2 )M 2 x (l-x + x 2 )M 2 1 



X X 



(B.25) 

and therefore the lower x-dependence don't affect the soft limits x — > 0, 1. 



B.3 The exact parametrization 



The Sudakov parametrization, as we will discuss in B.5 comparing our results with a full three body 
calculation, shows a bad behavior approaching x — > 0, namely when pi is soft. 

This because in (B.8) x cannot be considered exactly the variable describing the fraction of energy 
of the particle after the splitting process respect to its initial value. In order to correct this point, we 
need to introduce a different parametrization. Referring to the eikonal integral (B.5) we write: 



zE,k t ,0, 



;2£2 _ k 2 



xE,-k t ,0, x /x 2 E 2 -k* , . 



k 



Pi = 

p 2 = (yE,0,0-yE); 

and from energy and spatial momentum conservation we have: 

J x + y + z = 2, 

\ y/z 2 E 2 - k 2 t + v / x 2 E 2 - k 2 -yE = 0, 



together with the conditions: 

x 2 £ 2 > k 2 



t > 



z 2 E 2 >k 2 , <x,y,z <1. 



In this parametrization x can be considered exactly as the variable describing the fraction 
resolving the ambiguity noticed in the Sudakov parametrization. 

The scalar products appearing in (B.5) can be explicitly rewritten as: 



pi ■ k = xzE 2 + k 2 - ^z 2 E 2 - k 2 ^x 2 E 2 - k 2 , 
P2-k = 



Pi -P2 



yE 



yE 



zE + \ z 2 E 2 - kf 



xE+Jx 2 E 2 - k 2 



while for the phase space of the emitted soft particle we have: 

d 3 jfi dxdkl E 



(2vr)32p? 16vr 2 X 2 E 2 - kf 
In order to simplify the integration we note that: 



(B.26) 

(B.27) 
(B.28) 



(B.29) 

(B.30) 
of energy, 

(B.31) 
(B.32) 

(B.33) 
(B.34) 
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in the soft limit x — > from (B.29) it follows that z ~ y ~ 1, 



• since x 2 E 2 > A; 2 and < x < 1 it is possible to approximate xE 2 + k 2 ~ xi? 2 . 
As a consequence the scalar products in (B.31,B.32,B.33) are simplified: 



E 



E 



xE 



x 2 E 2 - k 2 



xE + Jx 2 E 2 - kf 



2E l 



(B.35) 

(B.36) 
(B.37) 



and the eikonal integral, considering (B.30) and introducing the mass Mw as a physical cutoff for the 
k 2 — > singularity, reduces to: 

"^ 2 / 4 dk 2 



D 



X-5-0 



9 



4E _ a 2 2 
that shows the behavior: 



In 



sx 
4M^ 



+ 21n 1 + 



4M 2 



2 

x 4M 2 , 



(B.38) 



(B.39) 



when sx 2 /4 ~ M^/, vanishing correctly when x = 2Mw / \fs. 

The symmetry propriety of the eikonal integral allow us to generalize this result for k soft just 
considering the substitution x — > z ~ 1 — x and therefore the kinematical end-points on the x variable 
for the exact parametrization are: 



2M W , 2M W 
— =- < x < 1 p^-. 



(B.40) 



As a conclusion we obtain the leading behavior of integrated splitting functions in correspondence of 
the two kinematical limit in the exact parametrization: 



x 1 



2M W 
2M W 



-Pexact ~ -L(x), 
X 



exact 



1 -x 



L(l-x), 



(B.41) 
(B.42) 



with L(x) given in eq. (18). 



B.4 The Collinear Approximation 

The eikonal approximation allows to highlight the singular behavior of the improved splitting functions 
in the soft regions x — > 0, 1, as shown e.g. in Eqs. (B.41, B.42). In order to extract the entire structure 
of the splitting functions and to show the factorization proprieties of our model independent approach, 
we need to go one step further, introducing the CA; we discuss now its main features, having in mind an 
illustrative explicit example. Following [17] we add to the SM Lagrangian a vector boson Z' with mass 
M 2> Mw belonging to an extra U(l)' gauge symmetry and singlet under SU(3)c <S> SU(2)l <8> U(1)y- 



2S 



In order to simplify our discussion, let us suppose that the Z 1 couples only with left electron and 
neutrino 

Ant = IlZ'^L, with: L = (u L , e L f. (B.43) 

At tree level we have two possible leptonic decay channel T ee = Yi{Z' — > e^e^) and T uv = Y^iZ' — > 
vlMl), related through total isospin conservation Y uu — Y ee = Tb- Considering the process Z' — > vl^l, 
labeling with P = (E, 0, 0, E) and P = (E, 0, 0, —E) the two back-to-back momenta of the two 
neutrinos (with E = M/2), and indicating with Y%{Z' — > vlMl) the corresponding two-body decay 
width for the amplitude squared we have at Born level: 

l-MBornl 2 = f L Tr (? 7 ^7 V ) e;(Q)e v (Q). (B.44) 

We calculate now the effect of adding one weak gauge boson emission, focusing on the three-body 
decay width T3 related to the process Z'(Q) — > vl{pi)vl(P2) Zrik) and using the CA. 

The key point of this approximation is the following: in the high energy regime M 3> M\y the 
leading contributions to the three-body decay are produced by the region of phase space where the 
emitted boson is collinear either to the final fermion or to the final antifermion, and in this region the 
three-body decay width is factorized with respect to the two-body one. 

Introducing: 

*r 3 - lM J (2^(Q — Pi — P2 — ^) (2 t 3 ) T 2 ;, ( ff 2 ;. (2 g,, . (B-45) 

it's possible to show in a simple way how the CA works both considering the factorization of the 
amplitude squared and of the phase space related to the final state. 

Considering for definiteness the case in which the weak gauge boson k is emitted along p\ direction, 
we depict in Fig. 6 the two Feynman diagrams involved in the computation of the amplitude 




Figure 6: Feynman diagrams involved in the calculation of the amplitude M.3 = M.a + -Mb of 
the decay process Z'(Q) — » vl(pi)W(P2) Z (k) . 



Concerning the factorization of the amplitude squared we have, working for simplicity in the 
massless limit: 

1 



iM A 
iMb 



9 



2cw 2pi • k 

_9 L_ 

2cw 2p2 ■ k 
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(B.46) 
(B.47) 



where T = i/L7^e (U (Q); therefore when the gauge boson is emitted along the p\ direction p\ • k — > 0, 
and the first diagrams diverges while the second one is finite; as a consequence squaring the amplitude 
M.3 = M.A + -M-B and using Feynman gauge divergences appear in the square |-Myi| 2 and in the inter- 
ference term MaM*b + M*^Mb- Using the Sudakov parametrization B.2 for the external momenta 
it's straightforward to verify that: 

now we use the only technical point of CA; writing: 

we simply observe that |.Ma| 2 diverges as k 2 — > 0; in the trace we can therefore approximate k ~ 
(1 — x)P and p2 ~ P: all the terms excluded by these approximations, in fact, softens the divergence 
and can be neglected. As a consequence we have: 

\M A \ 2 -^^xtl-xfTrhTr*?), (B.50) 

ZC W K t \ / 

and the remaining trace is exactly the one obtained in eq. (B.44) related to the process with no gauge 
boson emission, showing our first step towards the factorization of the amplitude squared. 
Considering the interference term A4a-M*b + M\Mb it's possible to use the same trick writing: 



2 z 

M A M* B -J^r^Tr (rTT*p) , (B.51) 

"w 



2r 2 k 2 



and as a result we obtain the factorization of the amplitude squared for the three-body decay in the 
CA respect to the two-body one: 



- W . (B . 2) 



Referring to the eikonal approximation in eq. (B.3), we see that using CA it's possible to factorize 
the amplitude only respect to its squared. 

We can apply the CA also considering the three-body phase space of T3 in eq. (B.45). Following 
[17] we therefore obtain a complete factorization: the three-body decay width can be expressed as the 
product of the two-body one times a collinear factor: 

(V , —7\~Ar (v> , -\ "2 1 l+x 2 dxdk 2 
dT 3 (Z -+ v L v L Z) « dT 2 (Z -> VLVL ) - — ^ — (B.53) 

Integrating over the final phase space we finally find: 
with s = M 2 . 
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As commonly done in parton models, we interpret the factor multiplying the two-body decay width 
as the parton distribution for finding a "neutrino parton" in the neutrino: 

where up to this point this expression take into account just the effects related to the real gauge boson 
emission. Adding a tree level delta term to account the process without electroweak emission and 
introducing virtual corrections at the same perturbative order of the real ones we obtain: 



where: 



2 s 2 



4(1 - s' w ) 



> (B-56) 



A complete set of integrated splitting functions is collected in Table 1. eq. (B.56) represents a 
concrete one-loop example, obtained through a direct calculation, of the logarithmical structure of the 
electroweak fragmentation function Dj^j. At this point one might be skeptical about the effective 
validity of the CA. In order to remove all doubt, we compare in B.5 our improved CA with the full 
result of a complete three-body calculation done in the context of the Minimal Dark Matter model [25], 
finding an excellent agreement. The one- loop fragmentation functions for the entire SM are discussed 
into Appendix C, in the more general context of the electroweak evolution equations [18]. 

At this point we are ready to evaluate in this example the energy spectrum of stable SM particles 
produced by the Z' decay and for definiteness we consider the case of the neutrino. We get: 

dN DM ^ UL _ 1 f dT 3 (Z' ^ u L VZZ) | dT 3 (Z' ^ u L e+W-) \ 



dx r to t i dx dx 

where r tot = T ee + T vv = 2Tb- From eq. (B.55) we have: 
dT 3 (Z' -> v L vZZ) 



dx 

and, in a similar way: 



T 2 (Z' u L VZ) D UL ^ L (x) = T B D UL ^ UL {x), (B.59) 



— = T 2 (Z ele L ) D £ -^ Ul (x) = F B D e -^ UL (x), (B.60) 

with D e -^ UL as in (C.2); finally we obtain the neutrino spectrum at perturbative order a 2 - 



(x)+D e -^ UL (x)]. (B.61) 



dx 2 

B.5 Full computation in the Minimal Dark Matter model 

In order to validate the eikonal and collinear approximations, we compare its results with a full compu- 
tation performed in a specific predictive model. We consider the weak corrections to DM annihilations 
into W + W~ as predicted by "Minimal Dark Matter" models, where DM only has electroweak inter- 
actions [25] . This generic situation is realized in the region of the MSSM parameter space where DM 
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could be the neutral component of the fermionic wino weak triplet with a value of the mass M ~ 3 TeV 
dictated by the cosmological DM relic density [25, 29]. The same situation can be realized with scalar 
DM and/or with DM lying in different representations of the weak group. Particularly interesting are 
two cases — a fermionic 5plet and a scalar 7plet with zero hypercharge — where DM is automati- 
cally stable because SM particles do not have the quantum numbers that could couple to such DM 
multiplets. Our computation applies in all such cases: at leading order we have the process DM DM 
->. W + W~ , and at NLO the three-body annihilation channels DM DM -)• W + W~Z and DM DM 
— > W + W~^ open up. The full expressions for the spectra for the emitted 7 and Z are quite lengthy, 
and therefore, we write them approximating Mjy ~ Mz and defining e = Myy/M <C 1, neglecting 
terms of C(e 2 ): 7 



dNpM- 
dx 



dNpM^z 
dx 



7T 
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(B.62a) 
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(B.62c) 



6(1 — x)x 



+ 



(1 — x)x 2 



ln(l - x) 



In the W spectra the first term is the leading order annihilation, the second one arises from 3-body 
processes with an additional Z, and the third term from processes with an additional 7. Eq. (B.62a) 
agrees with the result computed for fermionic DM in [30]. 8 

We now compare the full result with its collinear approximation, where the same quantities in 
Eqs. (B. 62a, B. 62b) are described through electroweak fragmentation functions. In the previous 



incidentally we find that at this order the result is the same for both scalar and fermionic DM. Notice that 
this Taylor expansion in e is only valid at x ^> e, and consequently, like the Sudakov approximation, does not 
correctly describe the kinematical boundaries. 

8 Note that the diagram there called "QED Internal Bremsstrahlung" in our language is ordinary EW 
bremsstrahlung from the initial DM and it does not give any log-enhanced effect because DM is non-relativistic. 
We also agree with [30] with the terms suppressed by AM/M\y, that we do not show because a AM would 
need a DM coupling to the Higgs, and consequently extra 2 body annihilations. 
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Z from W riidiiilion, M - 3 TeV 



y from W radiation 




x=E/M 




X 



Figure 7: Left, Z from W radiation: Comparison between our full result in the Minimal Dark 
Matter model (continuous yellow line), with its limit for e = Mw/M — > (blue dot dashed) 
and with our improved eikonal approximation (red dotted for the Sudakov parametrization and 
green dashed for the exact one). We show also the comparison with the naive standard partonic 
approximation (black continuous line). Right, 7 from W radiation: comparison between our 
full result (continuous red/blue line) with our improved splitting approximation in the exact 
parametrization (red/blue dashed) and the standard partonic one (red/blue dotted). 



sections we derived the results treating the phase space either within the Sudakov parametrization 
[Eqs. (B.16,B.17)] or exactly [Eqs. (B.41,B.42)]. We therefore compare the full 7 spectrum in 
eq. (B.62a) with: 



dx 



7T 



1 



-L(l-x) + 
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+ x(l — x) In 



M 2 wi 



and the full Z spectrum in eq. (B.62b) with: 
dN DM ^Z «2 



dx 



7T 



2c^ 



1 



L(l-x) + 



x 



L(x) + x(l — x) In 



M 2 W 



(B.63) 



(B.64) 



where L(x) is given in eq. (B.18) for the Sudakov parametrization and in eq. (18) for the exact phase 
space. Results are shown in Fig. 7, where we depict also the curve corresponding to the naive partonic 
approximation, where the upper bound on /j, 2 is chosen to be the typical hard scale of the problem s, 
obtaining a universal logarithmical factor lns/M^, that multiplies the splitting function. 

At small values of x the Sudakov parametrization (red dotted line in the left panel) shows a bad 
behavior compared with the full calculation and don't respect the correct kinematical end-points: 



2M W 



< x < 1 



2M W 



(B.65) 



This is because in the Sudakov parametrization, considering the splitting process i — > / + /', the 
variable x don't correspond exactly with the fraction of energy carried by the final particle / respect 
to its initial value. On the contrary the exact phase space (green dashed line) gives of course the 
correct kinematical boundary of the splitting process (B.40) and shows a correct agreement with the 
full calculation. 
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C One loop Electroweak Fragmentation Functions 



In the following we collect one loop EW fragmentation functions obtained solving the EW evolution 
equations in [18] for the entire Standard Model particle spectrum. 

EW evolution equations have been already constructed in [18], exploiting the SU(2)l symmetry, 
and classifying the states looking to their total isospin quantum numbers; consequently one have to 
apply a projection technique (explained in details in [18]) in order to convert them to their QCD-like 
formulation, i.e. labeling splittings with particle names: -Dj-s^ is the single leg fragmentation functions 
that encode the probability for a single initial particle i to become a final particle j. 

Since DM gives particle-antiparticle pairs, to reduce the combinatorics, we combine them into 
pairs of primary back-to-back SM particles I instead of a single particle i. A simple formula allows us 
to switch from the 'single leg' to the 'double leg' convention: 

Dj^j = C (Di-yj + I), q + I), . • + I), .,-) , (C.l) 

where c = 1/2 for complex final particles (such as the W or the v) while c = 1/4 for real ones (such 
as the Z or the 7). 



C.l Splitting of fermions 

We start considering DM that produces two back-to-back fermions, // with center-of-mass energy y/s, 
and compute the resulting partonic spectrum of other SM particles A, Df^A(z), with z = Ea/2\/s. 
As the // pair produces A and A with the same energy spectrum, we always average over particle 
and its anti-particle both in the initial and in the final state, even for real particles such as the Z. 
The initial fermion can be / = {e#, e^, ul, UL,d F , ur, (Ir\ and is identified by its T3 = { — 1/2, 0, 1/2}, 
its electric charge Q, its generation number. We define the usual coupling to the Z, gj = T3 — s^Q 



and to the photon, a em = s 2 v a2- We define the top Yukawa coupling yt 
and at = y^/Air. Neglecting all other fermion masses we get: 



mt/v with v ~ 174 GeV 



Dj^fix) = 5(1 
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D f ^z T (x) 
D f ^ WT (x) 
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- — f - P F ^ v (x), 
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2tt 
a 2 

2^ 



Q Pf^v(x), 
2T 2 P f ^ f (x), 
2T 2 P F ^ v (x), 



where /' is the SU(2)^ partner of / (e.g. /' = ei for / 
first line means that a fraction of the initial / with x ■ 
x < 1. 



+ (C.2a) 
(C.2b) 

(C.2c) 

(C.2d) 
(C.2e) 
(C.2f) 

vl and viceversa). The virtual term in the 
1 disappears into / or other particles with 
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All splittings of quarks are negligible with respect to the QCD splittings (not written) . For the t and 
b quarks there are extra splittings (to be summed to the ones listed above due to gauge interactions) 
due to the top quark Yukawa interaction, which also have minor effects. Again these splittings depend 
on the quark polarization; for simplicity, we write them for average unpolarized t, b quarks: 
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(C.3a) 
(C.3b) 
(C.3c) 
(C.3d) 



C.2 Splitting of Higgses 

The Higgs doublet H contains the physical Higgs h, as well as the goldstone components that describe 
the longitudinal polarizations Wl and Zl of the SM massive vectors. Splittings induced by the top 
quark Yukawa coupling have a significant effect; we describe them without specifying the polarizations 
of the t, b quarks, which make a negligible difference. 
For the physical Higgs h we have: 
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The same expressions hold for an initial Zl- For the longitudinal Wl we have: 
D Wl ^ Wl (x) = 6(1 -x 



3a 2 + ay ir 
J- "1 ~ : -if 



D 



(x) 



6(1 -x 
a 2 1 



2tt-4 



1 + — px* 



S^S 



+ — —Ps^s(x) + 



27T-4 



pvir,Yuk 



. a em j-. f \ 1 

-g^r p s^s(x) + 



2vr 4 



Ps-*s(x) = Dw L ->z L (x), 



Dw L ^z T (x) 

Dw L ->~f(x) 
Dw L ^ t (x) 



02 9eL 

2vr c 2 



Ps^y(x), 



ftp 



Ps->v(V), 



2tt 



3 at 



P, 



Yuk 



5->F 



(X). 



Eq. (C.5b) shows the QED corrections only, in case one needs to subtract them. 



C.4a) 

C.4b) 
C.4c) 
C.4d) 

C.5a) 
C.5b) 
C.5c) 
C.5d) 
C.5e) 
C.5f) 
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C.3 Splitting of vectors 

For the transverse W ± we find: 

D\y T ~^w T (x) = 6(1 — x 

Dw T ^z T (x) 
D Wt ^ 1 (x) 



1 i Q2 9 pvir 

2tt su ( 2 ) 



+ ^ Pv^v(x), 



^ Cw P V ^ V (x), 
Z7T 



°2 2 p / \ 

— s w P v ^ v (x), 



D 



D 



W T ^h 



(x) 
(x) 



(X2 1 

2tt 2 

a 2 1 

27 4 

a 2 1 

27 2 

and for the transverse Z we find: 



D Wt ^ Wl (x) 



N c P v ^ F (x), f = {e, i/ e , d, u; //, s, c; r, i/ T , b, t}, 

Pv^s(x) = D Wt ^ Zl (x), 

Pv^s(x), 



Dz T ^z T (x) 

Dz t ^W t ( x ) 
Dz T ^f(x) 

Dz T ^h(x) 

Dz T ->w L (x) 



6(1 - x) 



«2 



1 j. lo^ pvir i -W pvir 
1 + ^ ^ C W *SU(2) + "j" *U(1) 



«2 

27 

2^ 

2 

Py_, S (x) = D Zt ^Z l (x), 



2 C 2 ^ Py_>y(x), 



2 iV c Py^ F (a;), / = {v e ,e L ,e R ,u L ,u R ,d L ,d R , . . .}, 



2-kc 



w 



02 2g? 

2vr c 2 



Py^s(x), 



IV 



where we defined N c = 1 (3) if / is a lepton (a quark) doublet; N gen = 3, and: 



pvir , 1 pvir , at pvir 

- r su(2) — 2 v~* v 4 V ^ S gen y~>F-> 

P$i) = N gCQ (2Yi + Yi + 6Y^ + 3Y§ + 3Yl > )P^ F + 2YiP^ s . 

The 7 contributes to P\v T ^w T an d can be excluded by dropping 1 = c^ + s 2 ^ — > c 2 ^ in front of Py_>y, 
both real and virtual. (Notice that PYTHIA does not include QED radiation from W^). 
Finally for the photon we have: 



(C.6a) 
(C.6b) 
(C.6c) 
(C.6d) 
(C.6e) 
(C.6f) 

(C.7a) 
(C.7b) 

(C.7c) 
(C.7d) 
(C.7e) 



(C.8a) 
(C.8b) 



D 



7^7 
D-y^Wr (x) 

D ~t->W L (x) 
D^ f (x) 



(x) = 



6(1 - x 
17 



1 + 



o pvir , pvir 
2?r V 1 ^SU(2) + *U(1) 



2vr 
~2tT 



2 Py^y(x), 
2 Py^ S (:r), 

2Q 2 N c P v ^ F (x), f = {e L , e R , u L , u R , d L , d R , . . .}. 



(C.9a) 
(C.9b) 
(C.9c) 
(C.9d) 
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